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FILTERING  ANALYSIS  OF  A  DIRECT  NUMERICAL  SIMULATION 
OF  THE  TURBULENT  RAYLEIGH-BENARD  PROBLEM 


T.  M.  Eidion1 
M.  Y.  Hueaaini1 
T.  A.  Zang 

ABSTRACT 

A  filtering  analysis  of  a  turbulent  flow  hat  been  developed  which  provides  details  of 
the  path  of  the  kinetic  energy  of  the  flow  from  its  creation  via  thermal  production  to  its 
dissipation.  A  low-pass  spatial  filter  is  used  to  split  the  velocity  and  the  temperature  field 
into  a  filtered  component  (composed  mainly  of  scales  larger  than  a  specific  site,  nominally 
the  filter  width)  and  a  fluctuation  component  (scales  smaller  than  a  specific  size).  Variables 
derived  from  these  fields  can  fall  into  one  of  the  above  two  ranges  or  be  composed  of  a  mixture 
of  scales  dominated  by  scales  near  the  specific  size.  The  filter  is  used  to  split  the  kinetic 
energy  equation  into  three  equations  corresponding  to  the  three  scale  ranges  described  above. 

The  data  from  a  direct  simulation  of  the  Rayleigh-Benard  problem  for  conditions  where 
the  flow  is  turbulent  is  used  to  calculate  the  individual  terms  in  the  three  kinetic  energy 
equations.  This  is  done  for  a  range  of  filter  widths.  These  results  are  used  to  study  the 
spatial  location  and  the  scale  range  of  the  thermal  energy  production,  the  cascading  of 
kinetic  energy,  the  diffusion  of  kinetic  energy  and  the  energy  dissipation.  These  results 
are  used  also  to  evaluate  two  subgrid  models  typically  used  in  large-eddy  simulations  of 
turbulence.  Subgrid  models  attempt  to  model  the  energy  below  the  filter  width  that  is 
removed  by  a  low-pass  filter. 


1Research  was  supported  by  the  National  Aeronautics  and  Space  Administration  under  NASA  Contract 
No.  NAS1-18605  while  the  authors  were  in  residence  at  the  Institute  for  Computer  Applications  in  Science 
and  Engineering  (ICASE),  NASA  Langley  Research  Center,  Hampton,  VA  23665. 


1  Introduction 


The  numerical  simulation  of  turbulent  flows  has  provided  considerable  insight  into  flow 
physics  over  the  last  20  years[13,  14].  The  Rayleigh- Benard  problem  is  an  example  of  a 
specific  flow  which  is  better  understood  because  of  numerical  experiment*  [7].  The  high  cost 
of  direct  numerical  simulations  (DNS)  and  the  lack  of  confidence  in  extrapolating  subgrid 
models  in  large-eddy  simulations  (LES)  to  new  flows  has  prevented  this  technology  from 
being  more  widely  used. 

A  DNS  of  a  turbulent  Rayleigh- Benard  flow  was  computed  [8]  to  provide  a  database 
liar  a  more  extensive  analysis  of  a  turbulent  flow.  The  specifics  of  the  simulation  and  some 
preliminary  data  analysis  have  been  presented  previously [8].  The  focus  of  the  analysis  pre¬ 
sented  in  this  paper  is  to  use  low-bandpass  filtering  to  study  the  flow  field  generated  by  a 
direct  simulation.  The  intent  is  to  gain  insight  into  the  turbulent  flow  simulated  as  well  as 
to  understand  better  the  filtering  and  modeling  used  in  large-eddy  simulations. 

The  LES  approach  is  based  on  splitting  the  velocity  field  into  a  large-scale  (filtered 
velocity)  and  a  small-scale  (fluctuation  about  the  filtered  velocity)  component.  The  large- 
scale  component  is  calculated  from  the  appropriately  modified  Navier-Stokes  (and  thermal 
energy)  equations  (LES  equations)  and  is  resolvable.  The  modification  is  basically  done  by 
adding  a  subgrid  model  for  the  effects  of  the  small-scale  motion.  A  rigorous  derivation  of 
the  LES  equations  proceeds  by  spadally  filtering  each  term  in  the  original  equations  with  a 
low  bandpass  filter.  Details  of  the  derivation  of  the  LES  equations  can  be  found  in  [14].  The 
result  is  a  set  of  differential  equations  for  the  filtered  variables  (usually  the  velocities  and 
the  temperature)  with  any  terms  involving  a  fluctuation  variable  replaced  by  an  appropriate 
subgrid  model. 

To  gain  insight  into  the  LES  approach  the  total  dependent  variable  fields  can  be  computed 
using  the  DNS  approach  and  the  terms  involving  fluctuation  quantities  in  the  LES  equations 
calculated  to  evaluate  proposed  models.  This  has  been  previously  done  for  several  flows  for 
which  the  fluctuation  terms  (or  "turbulent  stresses”)  were  directly  calculated  and  compared 
to  models  by  a  correlation  analysis,  for  example  [1,  4].  In  the  present  study  a  correlation 
analysis  similar  to  those  of  previous  studies  is  done  with  the  following  differences.  Since  the 
vertical  direction  is  inhomogenous,  the  correlations  are  only  integrated  over  the  horizontal 
directions.  Thus  the  correlation  coefficients  are  functions  of  the  vertical  coordinate.  Since  the 
correlation  analysis  does  not  compare  the  mean  value  of  the  terms  of  interest,  the  horizontal 
average  of  these  terms  are  calculated  using  both  a  model  and  the  "exact”  turbulent  stress. 
These  horizontal  averages  are  compared  directly  to  evaluate  the  model[ll].  Also  the  fine 
grid  data  is  frequently  "smoothed”  onto  a  coarse  grid  with  a  grid  size  approximately  equal 
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to  but  smaller  than  the  filter  widths  used.  This  is  not  done  in  this  study  in  order  to  separate 
the  effect  of  the  numerical  method  from  the  effect  of  the  filter  width. 

While  insight  into  the  LES  technique  is  one  motivation  for  this  study,  the  filtering  analysis 
is  independent  of  LES.  A  range  of  filter  widths  from  the  grid  size  (or  viscous  scales,  whichever 
is  larger)  up  to  the  size  of  the  largest  flow  scales  of  the  problem  can  be  used  to  analyze  the 
data  from  a  particular  flow  simulation.  The  terms  in  the  various  kinetic  energy  equations 
based  on  the  filtered  and  fluctuation  velocities  can  be  calculated  as  functions  of  the  filter 
width.  This  yields  information  similar  to  that  of  energy  spectra  is  obtained  but  with  a 
different  perspective. 

The  filter  analysis  results  presented  here  are  calculated  from  the  output  of  a  DNS  of  a 
turbulent  Rayleigh- Benard  flow  (natural  convection)  at  a  Rayleigh  number  (Ra)  of  3.8  x  108 
and  a  Prandtl  number  (Pr)  of  0.76  .  These  values  were  chosen  because  of  the  availability  of 
experimental  and  other  simulation  data  at  similar  values.  Details  of  the  flow,  basic  variable 
defintions,  and  the  numerical  simulation  can  be  found  in  Appendix  A  and  in  a  previous 
publication  [8]  where  the  simulation  results  were  directly  compared  to  other  available  data, 
both  experimental  and  numerical[5, 2,  9]. 

Briefly,  the  simulation  entailed  the  computation  of  the  3-dimensional  velocity  and  the 
temperature  fields  on  a  128x64x64  grid.  The  ratio  of  horizontal  to  vertical  length  for  the 
computational  volume  was  4  and  2  in  the  xj  and  xa  directions,  respectively.  After  a  steady 
state  flow  was  developed,  data  was  collected  over  a  time  period  equal  to  10 /Wc,  where 
Wc  =  {NuPrRa)1^3  .  This  period,  which  should  consist  of  several  large  eddy  turn-over 
times,  was  found  adequate  by  Eidson[7].  The  average  Nusselt  number  (Nu)  calculated  from 
the  simulation  data  is  6.6  . 

The  analysis  reported  previously  [8]  shows  that  the  simulation  output  is  well  resolved. 
The  global  quantities,  the  Nusselt  number  and  the  RMS  averages  of  the  velocity  and  temper¬ 
ature,  are  in  good  agreement  with  the  experimental  results.  The  Nu  is  the  closest  of  several 
LES  and  DNS  studies  to  the  average  of  several  experimental  studies  {Nu  =  6.0).  The  vari¬ 
ation  in  x3  (the  vertical  coordinate)  of  the  various  velocity  and  temperature  quantities  and 
the  terms  of  the  total  kinetic  energy  equation  are  in  good  agreement  with  the  experiments. 
Most  encouraging  is  the  comparison  with  variance  and  skewness  measurements  of  the  tem¬ 
perature  and  its  vertical  derivative  by  Carroll[2].  While  the  agreement  in  magnitude  is  only 
fair,  the  changes  in  slope  with  x3  are  excellently  predicted.  These  results  along  with  some 
consistency  checks  for  the  volume  averaged  terms  in  the  various  kinetic  energy  equations 
which  are  made  in  the  current  analysis  provide  a  high  degree  of  confidence  in  the  simulation 
results. 
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2  Derivation  of  the  Filtered  Kinetic  Energy  Equa¬ 
tions 

The  choice  of  filter  to  separate  the  large-scale  and  small-scale  components  is  important.  The 
fundamental  structures  or  basis  functions  of  turbulence  are  not  sufficiently  understood  to 
select  the  best  filter  to  separate  large-scale  from  small-scale  flow  structures.  Sharp  cut-off 
filters  have  been  used  with  some  success  for  LES  studies,  but  they  do  not  allow  for  flow 
scales  near  the  filter  width  to  be  treated  separately  from  the  smaller  scales.  Consistent  with 
the  current  understanding  of  turbulence  any  smooth  low  bandpass  filter  is  satisfactory.  A 
Gaussian  filter  function  (shown  for  1  dimension)  is  generally  used  and  thus  was  selected  for 
this  study. 

ti(x)  =  f  G(z  —  x)u(x)dx  (1) 

J—OO 

G(i)  = 

7  =  6  (2) 

G  is  the  filter  function  and  A/  is  the  filter  width.  The  above  filter  is  applied  only  to  the  two 
horizontal  directions,  xj  and  x3,  to  avoid  the  complications  of  filtering  in  an  inhomogeneous 
direction,  x3.  The  fluctuating  component  (denoted  by  a  prime)  of  a  variable,  e.g.  it(x,t),  is 
defined  as  follows: 

u'(x,  t)  =  u(x,  t)  —  ti(x,  t )  (3) 

Five  filter  widths  are  used  to  analyze  the  simulation  data — 1/16, 1/8,  1/4,  1/2  and  3/4. 
For  comparison,  the  length  of  the  sides  of  the  control  volume  for  the  simulation,  Lm.,  are: 

K  =  4 
LXi  =  2 

=  1  (4) 

The  grid  dimensions  in  the  horizontal  directions  are  both  1/32  .  In  the  vertical  direction  a 
non-unform  Chebyshev  grid  is  used  with  the  average  grid  size  being  1/64  . 

The  derivation  of  the  filtered  kinetic  energy  equations  proceeds  as  follows: 

1.  The  time-dependent  partial  differential  equations  for  continuity  and  momentum  (Equa¬ 
tions  78  and  79  in  Appendix  A)  are  filtered  term  by  term. 
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2.  Filtered  velocity  product*  are  expanded  as  follows: 


Vi Vj  =  1^+vfa +*(%+*$ 

—  My + ~Lij + ~Cij + 32 $ 

(«) 

&ii  =  uiu'j 

(6) 

Cn  =  «&•  +  tiiuj 

(7) 

J 

III 

s' 

(8) 

II 

fcj 

C: 

1 

£ 

(») 

3.  Filtered  derivatives  of  a  variable  can  be  shown  to  be  equal  to  derivatives  of  the  filtered 
variable. 


~5u  da 

dx±  dx^ 

4.  Using  2.  and  3.  to  simplify  the  result  of  1.,  filtered  equations  are  formed. 


(10) 


(11) 


dMjj  _  _dP_,p 
di  dxj  dxi  T  dxj 

i+^±M  +  PrRaf6a  (12) 


5.  The  filtered  equations  are  subtracted  from  the  instantaneous  equations  to  form  fluctu¬ 
ation  equations. 


(13) 


-£+iv*i 

OXi  UXj 

+  *k  +  PrRaria 

uXj  OXj  OXj 


(14) 


4 


6.  A  kinetic  energy  equation  is  derived  from  a  dot  product  of  a  velocity  times  the  appro¬ 
priate  momentum  equation.  The  total  equation  is  the  total  velocity  times  the  instan¬ 
taneous  momentum  equation.  The  filtered  equation  comes  from  the  filtered  velocity 
and  the  filtered  momentum  equation  and  the  fluctuation  equation  uses  the  appropriate 
fluctuation  quantities.  The  cross  kinetic  energy  equation  is  formed  as  the  fluctuation 
velocity  times  the  filtered  momentum  equation  added  to  the  opposite  (filter  velocity 
and  fluctuation  momentum  equation). 

•  Total  Kinetic  Energy  Equation 

SB  _  8Mg  +  P)]  .  -  dfujSjj) 

at  ~  ezj  *  dxj 

+PrRa(utT)-lprS,iSij 


=  D’  +  Pr-t 

(18) 

Bs  r** 

(16) 

, - 9{*j(E  +  P))  „_«M«) 

D  -  e*i  +Pr  Hi 

(17) 

PT  =  PrRafaT) 

(18) 

6  =  \PrSi’Sii 

(19) 

•  Filtered  Kinetic  Energy  Equation 
dEF  _  dfajEr  +  P)} 


9. 


j)T  2  jy*r  _  +  %f  +  Thj) 

dxj 

W'sWff  +  Wg  +  Wl' 

p*  =  -ptf-pg-p*f 

if  =  PrRa'faT) 

f  =  \prSij5ii 

•  Crow  Kinetic  Energy  Equation 

9Bc  _  d^P  +  W)  .  n_a(u'?„-Hi.S!,) 

9t  dxj  dxj 

.  -  ~  fij)  __  . 

™  dXj  ‘  0*, 

■  ~  •%>)  _  /^Ssi 

ax,  *  a*,- 

- 1 {—&  +  PrRa'^T  +  SjT') 
—PrSijS'ij 

=  £?,c7  +  ^  +  if  -  6C 
=  Ec  +  Pc  +  P?-ec 

Ec  =  u-tii 

q*c g  W£±*£) + pr^u^+^ 

5x,‘  dx^ 

Dc  =  +  appropriate  stress  diffusion  terms 


+«  -  «'/') 


Pc  =  -(P‘c-Pl!,)-Pil-(PS-Pkl)-Pfl 
-{PSl  -  Pi ') 
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if  =  i>rJ!a'(uy'  +  SsT') 

(34) 

(35) 

•  Fluctuation  Kinetic  Energy  Equation 

m  ,  P,  A»X) 

at  dxj  T  Bx, 

1  v!  ~  £*il  I  ~  ^ij) 

‘  0x,  +u<  dxj 

FST 

H^  +  PrRaXvp1) 

T*Wr 

=  IT* +  W* +  />*-£* 


=  P5  +  P,  +  P/-es 

(36) 

B,  =  i«X 

(37) 

“  dxj  1  a*, 

(38) 

2?5  =  D*s  +  appropriate  stress  diffusion  terms 

(39) 

Ws  =  (W*-W#)  +  (W*-Wg)-Wl' 

(40) 

p*s-{p*-p!.<)-(p*-pX)  +  pI' 

(41) 

Pi  =  +PrRa'(x^r) 

(42) 

<?S\  PrSSM 

(43) 

The  production  terms  are  composed  of  the  thermal  production, 
production,  Pf.  The  velocity  production  takes  the  following  form: 

Pt  and  the  velocity 

■pah  _  a  a  rrb 

*Z  -  "ij  ) 

(44) 
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where  Zy  ia  one  of  the  "turbulent  stresses”  —  the  Reynolds  stress,  Rij}  the  Gross  stress,  C,,, 
or  the  Leonard  stress,  Zg.  Another  term  which  appears  in  the  analysis  with  a  role  similar 
to  a  stress  is  Af,-,-.  The  superscript  b  refers  to  the  filtered  value  of  the  stress(6  =  /)  or  to  the 
unfiltered  or  total  stress(&  =  t).  These  terms  are  called  production  terms  since  they  have 
the  role  of  extracting  energy  from  the  flow  at  one  velocity  scale  and  adding  energy  to  the 
flow  (or  producing  energy)  at  another  scale.  The  term  5^°  is  the  strain  tensor  which  can  be 
formed  from  filtered  velocities  (a  =  /),  or  from  fluctuation  (or  perturbation)  velocities 
(a  =  p),  S'ij.  The  notation  for  each  of  the  particular  production  terms  that  appear  in  the  set 
of  filtered  kinetic  energy  equations  derived  in  this  paper  is  given  in  Appendix  B. 

The  production  term*  are  only  one  form  that  the  work  of  stresses  can  take  in  a  kinetic 
energy  equation.  Another  form  is  as  follows: 


(4!) 


This  form  will  simply  be  called  the  work  term.  The  work  term  can  be  expanded  into  the 
production  term  plus  another  work  term,  Tjf  which  measures  the  diffusion  of  the  work  of 
the  stresses. 


wt  = 

*Tf-P? 

(46) 

Tf  = 

aK’Zt) 

dXi 

(47) 

Tj *  will  be  referred  to  as  the  turbulent  stress  diffusion.  The  volume  average  of  the  production 
term  is  the  negative  of  the  volume  average  of  the  work  term  for  the  same  Z  and  (<*6).  Using 
the  straight-forward  derivations  for  the  kinetic  energy  equations  described  above,  the  work 
form  naturally  appears.  Equation  46  can  be  used  to  modify  pointwise  the  kinetic  energy 
equations  to  include  the  production  form.  In  general,  the  production  form  and  the  work 
form  are  not  equal  locally.  One  exception  is  the  following: 

Pg  =  WS  (48) 

Each  of  the  kinetic  energy  equations  can  be  written  as  "the  time  derivative  of  a  kinetic 
energy  equals  diffusion  terms  plus  production  (or  work)  terms  minus  a  dissipation  term.” 
If  the  flow  has  reached  a  "turbulent  steady  state”  then  the  time  derivative  of  the  kinetic 
energy  averaged  over  a  suitably  large  spacial  region  should  be  zero.  The  diffusion  terms  are 
of  the  form — the  sum  of  the  spacial  derivative  in  each  of  the  3  directions  of  a  group  of  terms. 
The  volume  average  of  the  diffusion  terms  can  be  shown  to  equal  zero  for  the  boundary 
conditions  in  this  simulation;  specifically,  periodic  boundary  conditions  in  the  horizontal 
directions  and  zero  velocities  at  the  walls.  Otherwise,  these  terms  are  not  zero  even  for  a 
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horizontal  average.  The  disappearance  of  the  time  derivative  and  diffusion  terms  for  suitable 
averages  partially  explains  why  a  "production  equals  dissipation”  assumption  is  often  made 
in  turbulent  theories[12, 10].  All  the  terms  in  these  kinetic  energy  equation  are  evaluated  in 
this  study. 

The  filtered  momentum  and  continuity  equations  become  the  LES  equations  after  the 
appropriate  terms  involving  fluctuation  variables  have  been  modeled.  Actually  the  fluctua¬ 
tion  kinetic  energy  equations  (and  even  equations  for  Rij  and  C#)  with  appropriate  models 
can  be  incorporated  into  the  LES  equations,  but  this  has  not  been  generally  done[6].  For 
the  present  study  Equation  15  will  be  called  the  total  equation  since  the  principal  dependent 
variable,  E,  is  a  function  of  total  velocities,  u<.  The  filtered  kinetic  energy  equation  also  can 
be  called  the  resolvable  scale  equation  since  the  principle  variable,  Ep,  is  a  function  of  the 
filtered  velocities,  which  would  be  directly  resolved  in  the  LES  technique.  The  perturbation 
kinetic  energy  equation  can  be  viewed  as  a  subgrid  equation  since  the  principle  variable,  Es , 
is  a  function  of  velocities  with  scales  below  the  filter  width  and  for  LES,  the  filter  width  is  of 
the  same  order  as  the  grid  size.  The  cross  kinetic  energy  equation  can  be  called  the  transfer 
scale  equation  since  it  has  been  suggested  that  Ec  is  principally  composed  of  velocity  scales 
that  transfer  energy  from  the  filtered  region  to  the  fluctuation  region[l].  Equations  20,  28 
and  36  will  be  referred  to  as  the  set  of  filtered  kinetic  energy  equations. 

The  total  potential  energy  equation  is  included  for  reference.  The  potential  energy  equa¬ 
tion  is  derived  by  multiplying  the  thermal  energy  equation  by  -PrRa(x3).  Note  that  the 
term,  PrRa(u3Tr)t  removes  energy  (has  a  negative  sign)  from  the  potential  energy  equation 
where  it  adds  energy  (has  a  positive  sign)  to  the  total  kinetic  energy  equation. 

¥ + =  -Pra,(“’r')  ~  PrRa  («) 

$  is  the  non-dimensional  potential  energy.3 

$  =  - PrRaTrx3  (50) 

aIn  dimensional  form,  4  is  defined  ss  follows: 

4  =  -  pe-ri  =  -gpTr*  j 

Pa 
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3  Data  Analysis  Details 


When  the  kinetic  energy  equations  are  horizontally  averaged,  the  cross  kinetic  energy  equa¬ 
tion  can  be  simplified.  The  simplification  is  a  consequence  of  the  following  identity. 

<fS>  =  y^/(*)(y^G(*-*')y(*')d®^dz 

=  jT  $(*')  (J_"  G(x  -  x')f(x)dx  j  dx' 

=  <}g>,  °°  (51) 


where  f  and  g  are  smooth,  continuous  functions  of  x.  The  <  >  denote  a  horizontal  average. 
When  the  above  identity  is  applied  to  certain  production  and  work  terms  the  following 
relationships  can  be  derived. 


OXj 


=  <Ui- 


SMt 


-Ui- 


dM, 


•; 


dx j  dlj 


_  SKI  a  SXXij 

OXj  OXj 

,SJXij 


-M>  = 

UXj  (JXj  OXj 

dUij  _  SEij 

,«L 

“  <li‘  dXj  > 

<  >=<  lOOi  > 

*  dXj  ’  dXj 


<  S&j .  ~ —  —  <  S-jMij  > 

<  >=<  S'ijUij  > 

<  -  Etf)  >=<  > 


(52) 


(53) 

(54) 

(55) 

(56) 

(57) 


These  relationships  can  be  used  to  simplify  the  horizontally  averaged  cross  kinetic  energy 
equation  to  the  following  form: 

<  TT  >=<  D*c  +  2Wc  +  2WR  +  2Wm  +Pt~*c>  (58) 

ot 
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Although  these  equations  are  derived  for  continuous  arithmetic,  the  derivation  for  discrete 
arithmetic  is  straight-forward. 

The  data  analysis  proceeded  as  follows: 

1.  The  velocity  and  temperature  values  at  each  grid  point  were  saved  for  322  time  steps 
during  a  Steady  state”  portion  of  the  DNS.  These  time  steps  are  equally  spaced  over 
a  time  period  of  \0jWc.  For  each  filter  width,  the  various  production,  work,  diffusion 
and  dissipation  terms  are  calculated  at  each  grid  point  for  a  particular  time  step.  All 
derivatives  and  multiplications  of  variables  are  done  in  the  same  manner  as  the  original 
simulation[8].  The  analysis  arithmetic  is  done  using  32-bits  only  compared  to  64-bits 
used  in  the  original  simulation. 

2.  The  terms  are  then  horizontally  averaged,  denoted  by  <  >,  and  volume  averaged, 
denote  by  {  }.  The  data  was  then  “long”  time  averaged  over  about  80%  of  the  "steady 
state”  portion  of  the  simulation — the  time  period  shown  in  Figure  1. 

In  the  plots  of  the  horizontal  averages  versus  x3,  the  curves  are  not  always  symmetrical 
about  x3  »  0,  particularly  the  production/work  terms  (e.g.  Figure  16).  Specifically,  some  of 
the  variations  with  x3  (or  "humps”)  appear  to  be  too  irregular  to  be  a  “general  steady  state” 
solution  that  would  be  found  regardless  of  initial  conditions  (and  possibly  small  changes  in 
boundary  conditions  such  as  the  aspect  ratio  of  the  computational  volume).  The  horizontal 
and  volume  averages  of  the  velocities  should  tend  to  zero  since  viscosity  should  dissipate  any 
mean  horizontal  motion.  "Non-symmetric  or  non-zero”  initial  conditions  or  a  weak  organized 
flow  caused  by  an  instability  may  persist  for  a  long  period  of  time[3].  Plots  of  the  horizonal 
velocity  averages  showed  that  a  small  patterned  flow  did  exist  during  the  entire  simulation. 
In  Figure  la,  {uj}  and  {ttj}  are  small  compared  to  the  RMS  levels  in  Figure  lb  but  show  no 
tendency  to  dissipate.  A  change  in  the  pattern  of  {u<}  was  observed  in  Figure  lb  between 
the  first  and  second  half  of  the  simulation.  Differences  in  <  u,-  >  versus  x3  between  the  two 
halves  are  also  observed.  The  pattern  in  Figure  2a  for  4  selected  time  steps  is  observed  during 
the  entire  simulation  with  the  peaks  oscillating  several  times  from  positive  to  negative.  This 
pattern  could  have  resulted  from  a  large  cellular  motion  that  filled  the  entire  vertical  region. 
The  pattern  in  Figure  2b  is  observed  from  the  middle  to  the  end  of  the  simulation.  This 
form  could  result  from  cellular  motion  of  about  half  the  size  of  the  distance  between  the 
plates.  <  u3  >  is  several  orders  of  magnitude  smaller  than  <  u3  >  and  <  u3  >  due  to 


the  effect  of  the  continuity  equation  end  the  periodic  boundary  conditions  in  the  horisontal 
directions.  These  facts  lead  one  to  suspect  that  the  component  of  the  flow  causing  the  non- 
•ero  magnitude  of  <  u\  >  and  <  uj  >  is  related  to  the  irregular  "humps'  in  the  velocity 
production  terms.  <iij  >  and  <  tij  >  are  a  maxi  mum  of  about  20  at  each  time  step.  The 
production  terms  are  proportional  to  velocity  cubed  and  the  "humps”  are  generally  of  order 
20s  =  8000.  Also  the  xj  sise  of  the  “humps"  is  consistent  with  the  variations  in  zj  seen  in 
Figure  2b  as  well  as  similar  plots  at  other  time  steps  (not  shown). 

One  other  source  of  inaccuracy  that  is  investigated  is  the  effect  of  aliasing.  Although 
there  is  no  apparent  aliasing  problems  in  the  simulation,  several  of  the  operations  in  this 
analysis — multiplication,  taking  derivatives,  extracting  perturbation  variables — increase  the 
relative  sise  of  the  higher  wavenumber  components  of  the  field  being  calculated.  As  a  check, 
the  original  velocity  fields  were  expanded  by  spectral  methods  onto  a  grid  with  twice  the 
resolution  and  a  limited  set  of  results  were  calculated  on  this  grid.  Although  some  differences 
in  the  results  are  observed  locally  (a  small  3-dimensional  spacal  region  at  one  time  step)  for 
the  most  sensitive  quanties  bong  calculated,  no  significant  difference  is  found  between  the 
horisontal  averages  calculated  on  the  original  grid  and  on  the  expanded  grid. 

4  Volume  Averages 

The  data  calculated  for  the  five  filter  widths  specified  previously  are  discussed  below  where 
the  terms  of  the  kinetic  energy  equations  have  been  averaged  over  the  computation  control 
volume.  The  filter  widths,  listed  in  Table  1,  can  be  compared  to  the  size  of  the  control 
volume  for  reference,  which  is  l(height)  by  4  by  2  in  units  of  non-dimensional  length.  The 
data  is  long  time  averaged. 

In  Table  1  the  non-zero  volume  averages  for  the  terms  in  Equations  20,  28  and  36  along 
with  the  kinetic  energy  associated  with  each  equation  are  shown.  Also  the  values  for  the 
diffusion  terms,  which  should  give  zero,  are  included  as  a  measure  of  the  numerical  accuracy 
of  the  original  simulation  and  current  analysis.  In  Table  2  the  fraction  of  kinetic  energy  in 
the  scales  below  the  filter  width  is  given. 

These  volume  averages  are  plotted  versus  filter  width  in  Figures  3a-d.  For  the  small  filter 
widths  the  filtered  kinetic  energy  (and  terms  in  the  associated  equation  )  dominates  with  the 
cross  kinetic  energy  (and  terms)  being  larger  than  the  fluctuation  equivalent.  At  the  larger 
filter  widths,  the  energy  is  more  evenly  divided.  In  each  of  the  3  kinetic  energy  equations 
the  thermal  production  and  the  dissipation  terms  dominate  and  almost  balance  each  other. 
The  other  non-zero  volume  averages,  the  production/work  terms,  are  smaller. 

The  production  of  the  fluctuating  flow  is  assumed  generally  to  be  due  to  thermal  effects 
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OOW) 

3264 

1639152 

1630896 

0 

filter 
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E, 

D*r 

p; 

e* 

Pr 

1/16 

8712 

3216 

1605264 

1559800 

•37386 

1/8 

8351 

3088 

1520776 

1403512 

•109760 

1/4 

7326 

2704 

1302627 

1080189 

-216349 

1/2 

5261 

1936 

916234 

657181 

•254166 

3/4 

3629 

1280 

632675 

418080 

-209019 

Ec 

Pt 

e* 

P° 

1/16 

135 

48 

32961 

67605 

34498 

1/8 

466 

160 

107777 

193472 

85342 

1/4 

1254 

480 

260544 

364362 

102924 

1/2 

2333 

912 

424360 

429592 

2853 

3/4 

2819 

1072 

488818 

406334 

-86228 

Es 

D*s 

PS 

e5 

pS 

1/16 

3 

0 

929 

3494 

2540 

1/8 

33 

16 

10600 

33913 

23207 

1/4 

270 

80 

75983 

186347 

110073 

1/2 

1255 

448 

298560 

544125 

244862 

3/4 

2401 

896 

517659 

806482 

287344 

Table  1:  Volume  Averaged  Kinetic  Energy  Equation  Terms 


filter  width 

( Eg+Es ) 

( Ep+Ec+Es ) 

1/16 

0.015 

1/8 

0.056 

1/4 

0.172 

1/2 

0.405 

3/4 

0.590 

Table  2:  Fraction  of  kinetic  energy  below  the  filter  width. 
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foe  the  large-scale  motion  and  due  to  the  velocity  rereading  mechanism  for  the  generation 
of  intermediate  and  small- scale  motion.  The  results  of  the  current  simulation  axe  able  to 
quantify  this  notion  at  least  for  one  Ra.  The  thermal  production  occurs  at  scales  larger  than 
1/8  based  on  the  observation  that  P$  +  fj?,  the  thermal  production  below  the  filter  width, 
is  reduced  to  7%  of  the  total  dissipation  at  this  filter  width.  At  the  largest  filter  width,  the 
equivalent  value  is  60%.  Thus  the  thermal  production,  while  concentrated  at  the  largest 
scales,  is  spread  over  a  significant  range. 

The  velocity  production  or  cascading  energy  is  shown  in  Figure  3c.  As  the  filter  width 
increases,  the  velocity  production  terms  in  each  of  the  3  kinetic  energy  equations  must 
approach  sero.  Assuming  that  these  terms  have  reached  their  maximum  value  near  the 
largest  filter  used  in  this  study,  only  about  20%  of  the  total  dissipation  (based  on  Ps  at 
A /  =  3/4)  is  involved  in  the  cascading  mechanism. 

In  Figures  4  and  5  simple  schematics  of  the  energy  flows  are  presented.  The  conversion 
of  energy  from  one  type  to  another  occurs  over  a  range  of  flow  scales  and  is  not  generally 
local  to  a  specific  special  region.  (The  effects  of  a  core  and  boundary-layer  region  are  not 
distinguished  in  this  volume  average  analysis.)  The  average  transfer  of  energy  from  the 
filtered  to  the  cross  to  the  fluctuation  equation  is  done  by  the  velocity  production/work 
terms.  For  the  flow  simulated  in  this  study,  energy  on  the  average  is  extracted  from  the 
filtered  equation  and  added  to  the  fluctuation  equation.  The  cross  velocity  production  was 
either  positive  or  negative  depending  on  the  filter  width. 

5  Horizontal  Averages 

In  this  section  the  variation  with  ®3  of  the  horizontal  averages  of  the  terms  of  the  set  of 
filtered  kinetic  energy  equations  (Equations  20,  28  and  36)  are  examined.  Plots  of  the  terms 
in  the  total  kinetic  energy  equation  (Equation  15  and  Figure  6a)  and  the  total  potential 
energy  equation  (Equation  49  and  Figure  6b)  are  included  for  reference. 

The  first  observation  is  the  existence  of  boundary-layer  regions  near  the  wall  where  the 
various  terms  have  significantly  different  behavior  from  the  center  layer.  In  a  previous 
analysis  of  this  simulation  data,  a  boundary-layer  region  of  thickness  0.23  next  to  each  wall 
was  determined[8].  This  agreed  with  experimental  data,  in  particular  that  of  Carroll[2].  All 
of  the  terms  studied  in  this  analysis  exhibited  a  region,  significantly  different  in  character 
from  the  core  region,  near  the  walls  with  approximately  this  thickness.  (Figures  7  to  9) 
In  the  core  region  most  of  the  kinetic  energy  equation  terms  are  relatively  constant  with 
the  major  exceptions  being  the  production/work  terms.  The  production/work  terms  are 
generally  smaller  in  magnitude  than  the  other  terms  in  the  same  equation  with  the  same  filter 
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width.  Their  variations  with  Xj  could  be  due  to  the  low  frequency,  patterned  Sow  mentioned 
previously  that  persisted  throughout  the  simulation.  In  the  boundary  layer  almost  all  the 
terms  are  significant  and  exhibit  large  variations  in  x3.  The  production  terms  approach 
aero  as  they  near  the  walls  while  the  dissipation  and  diffusion  terms  (due  to  the  viscous 
component,  see  below)  approach  a  maximum  value  at  the  wall. 

The  total  diffusion,  D,  in  each  of  the  three  equations  has  a  similar  variation  with  x3 
(Figure  7).  The  diffusion  is  transferring  energy  from  the  core  to  the  boundary  layer.  The 
diffusion  in  the  filtered  equation  increases  as  A/  is  decreased.  However,  as  A/  is  decreased, 
the  diffusion  in  the  fluctuation  equation  becomes  smaller  until  it  is  negligible.  The  diffusion 
in  the  cross  equation  is  larger  than  that  in  the  fluctuating  equation  but  also  decreases  with 
decreasing  A/.  The  diffusion  is  composed  of  the  stress  diffusion  plus  the  two  component  of 
D* — the  pressure  component  (pressure  plus  kinetic  energy  in  the  filtered  equation)  and  the 
viscous  component.  Both  contribute  mainly  near  the  wall,  although  the  pressure  component 
has  a  significant  core  contribution  in  some  cases.  Figures  8a  and  8b  are  examples  of  the 
typical  variation  with  X3  of  these  terms,  which  are  similar  to  the  same  terms  in  the  total 
kinetic  energy  equation  (Figure  6a).  The  stress  diffusion  (Figure  8c)  has  a  shape  similar  to 
the  pressure  component  and  is  generally  the  smallest  of  the  diffusion  components. 

An  example  of  the  dissipation  terms  is  shown  in  Figure  9a.  There  was  little  difference 
in  the  shape  of  the  dissipation  versus  x3  for  the  different  equations  and  filter  widths.  The 
thermal  production  terms  also  exhibited  this  consistency  of  shape.  An  example  is  shown  in 
Figure  9b.  The  magnitudes  of  these  terms  for  different  A /  can  be  inferred  from  Figures  3b 
and  3a.  One  other  quantity  of  interest  is  the  kinetic  energy  itself.  The  variation  with  x3  of 
the  three  kinetic  energies  was  similar  although  the  fluctuation  kinetic  energy  did  not  exhibit 
the  peak  in  magnitude  in  the  boundary  layer  as  did  the  other  two  (Figures  9c  and  3d). 
All  three  showed  no  significant  variation  in  shape  with  A/.  The  variation  with  time  of  the 
kinetic  energies  can  be  inferred  from  Figure  lb. 

6  Kinetic  Energy  Equations 

All  the  terms  of  the  filtered  kinetic  energy  equation  are  plotted  together  in  Figure  10  for 
3  filter  widths.  The  dissipation,  e,  is  positive  when  it  extracts  energy,  which  is  opposite  to 
the  other  terms — due  to  writing  —  e  for  the  dissipation  in  the  kinetic  energy  equations.  The 
terms  have  similar  variations  in  x3  to  those  of  the  terms  in  the  total  kinetic  energy  equation 
(Figure  6a), except  that  the  velocity  production  terms  are  non-zero  although  they  are  small. 
The  volume  averages  for  these  same  terms  except  for  the  diffusion  which  is  zero  are  plotted 
versus  A/  in  Figures  11a  and  lib.  While  the  thermal  production  balances  the  dissipation 
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globally  (volume  averaged),  this  is  not  true  locally  (horiaontally  averaged),  particularly  in 
the  core.  For  small  A /  the  dissipation  in  the  core  increases  but  still  does  not  balance  the 
production. 

For  a  large  A/,  the  filtered  equation  describes  the  energy  flow  of  the  large  scales.  r  .om 
Figures  12a  and  12b  the  energy  extracted  from  the  large-scale  flow  by  the  production  terms 
is  mainly  in  the  boundary  layer.  For  the  smaller  A /  in  Figure  12c  the  fluctuations  due  to  the 
scales  smaller  than  the  filter  width  are  extracting  energy  more  uniformly  across  the  entire 
fluid  layer.  Another  observation  is  the  change  in  importance  of  the  three  production  terms 
shown  in  Figure  lib.  At  small  A/,  P</  and  are  dominant.  At  the  larger  A/  where  C,-,- 
and  Zjj  approach  zero,  Pr  becomes  largest — although  this  term  also  must  approach  zero 
for  very  large  A/. 

The  cross  kinetic  energy  equation,  Equation  28,  contains  9  production/work  terms.  For 
horizontal  and  volume  averages,  a  simplified  form  reduces  the  number  of  terms  to  3  (Equa¬ 
tion  59).  Since  no  particular  physical  significance  can  be  associated  with  either  form,  the 
data  will  be  plotted  for  the  reduced  form  for  simplicity. 

The  cross  kinetic  energy  equation  terms  are  shown  in  Figures  13  and  14.  For  larger  A/, 
the  variation  of  the  terms  with  Z3  is  similar  to  that  in  the  filtered  equation  (Figure  1C). 
One  difference  is  that  the  dissipation  is  more  significant  in  the  core  for  the  cross  equation 
than  for  the  filtered  equation.  For  smaller  A/,  the  velocity  production  is  almost  equal  to  the 
thermal  production;  whereas,  in  the  filtered  equation  the  thermal  production  dominates  the 
velocity  production  for  all  A/.  This  is  due  to  the  large  magnitude  of  the  velocity  production 
component,  Pj/. 

Similar  plots  for  the  terms  in  the  fluctuation  kinetic  energy  equation  are  shown  in  Fig¬ 
ures  15  to  17.  To  simplify  the  plots,  the  two  velocity  production  terms  involving  Ca  have 
been  added  together.  The  same  was  done  with  the  Rij  terms.  In  each  case,  one  of  the  two 
terms  is  significantly  smaller. 
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The  data  for  the  fluctuation  equation  is  similar  to  that  for  the  cross  and  filtered  equation. 
However,  the  boundary  layer  and  core  regions  are  not  as  well  defined  as  in  the  other  two 
equations.  On  close  inspection  there  appears  to  be  a  region  in  the  core  merging  with  the 
boundary  layer —  0.2  <  |®3|  <  0.4  — where  the  velocity  production  and  the  diffusion  have 
increased  importance.  Another  difference  from  the  other  two  equations  is  that  the  thermal 
production  does  not  approximately  balance  dissipation  even  on  a  volume  averaged  basis. 
The  velocity  production  is  even  slightly  larger  than  the  thermal  production  for  small  A/. 
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The  magnitude  of  P*  is  not  mainly  due  to  one  component,  a a  for  the  cross  equation,  but 
the  “L,  C  and  R"  components  vary  in  importance  as  a  function  of  A/  (Figure  17b). 

7  Subgrid  Scale  Modeling 

Calculation  of  the  terms  in  the  kinetic  energy  equation  using  the  direct  simulation  data  can 
be  used  to  study  subgrid  models  for  LES.  For  example,  the  production  terms  for  Hij+Hij  can 
be  calculated  directly  (refered  to  as  “exact”)  and  compared  to  the  same  term  calculated  with 
a  model  for  22$  +7$.  The  "exact”  and  modeled  production  terms  can  then  be  compared 
using  a  correlation  analysis [1,  4].  However,  a  direct  comparison  of  some  average  of  the 
production  terms  (such  as  an  average  over  homogeneous  directions)  can  be  used  also  to 
evaluate  the  models.  Since  the  selection  of  models  and  model  constants  is  often  based  on  the 
correct  prediction  of  energy  transfer  between  different  flow  scales  in  many  LES  studies  [1],  the 
correct  prediction  of  the  magnitudes  of  the  production  terms  should  be  an  important  model 
criterion.  In  particular,  P^.c  =  Pr  +  Pc  ,  which  extracts  energy  from  the  filtered  flow 
field  (Equation  20),  should  provide  the  most  important  criterion.  The  analysis  in  this  study 
is  only  the  first  step  in  evaluating  subgrid  models.  Conclusions  drawn  from  this  analysis  of 
DNS  data  which  suggest  changes  in  subgrid  modeling  must  be  verified  in  an  actual  LES. 

The  main  thrust  of  this  study  was  to  evaluate  the  behavior  of  commonly  used  subgrid 
models  for  natural  convection  flows.  To  the  authors  knowledge,  most  studies  of  subgrid 
modeling  using  DNS  data  have  been  done  using  isothermal  flows.  The  model  comparisons 
were  made  using  a  range  of  filter  widths  from  0.0625  to  0.75  .  These  filter  widths  were  chosen 
so  that  the  split  of  kinetic  energy  between  the  filtered  velodty  field  and  the  subgrid  velocity 
field  would  include  the  range  found  in  typical  LES  studies.  The  Smagorinsky  model  and  the 
scale-similarity  model  were  the  models  tested  since  these  are  the  two  most  frequently  used 
models.  The  Smagorinsky  model  is  [7]: 

>55  -  (62) 

if  -  Wi 

•n 

S  = 

A s  =  minimum  (A j  or  distance  to  nearest  wall) 

Cr  —  model  constant  (63) 

The  scale-similarity  model  is  [1]: 
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filter  width  Ck 


1/16  0.047 

1/8  0.055 

1/4  0.069 

1/2  0.114 

3/4  0.175 


Table  3:  Smagorinsky  model  constant,  Ck,  adjusted  so  that  the  model  predicts  the  “exact" 
value  of  Pr+c  *t  each  filter  width. 


-  ^SijUkK  =  («i«i  -  «i«i)  -  ^u(tinUn  -  5»«n) 
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The  value  of  Ck  for  the  Smagorinsky  model  is  one  of  the  results  of  the  current  analysis. 
Generally,  Ck  =  0.21  has  been  used  in  natural  convection  simulations  [7],  but  this  value 
assumes  that  the  grid  size,  A,  is  used  as  the  length  scale,  A,.  Studies  for  homogeneous 
turbulence  have  shown  that  Ck  should  be  reduced  by  a  factor  of  2  or  3  when  A /  is  used 
instead.  Therefore,  a  value  between  0.07  to  0.14  should  be  assumed  as  a  target  value  of 
Ck-  If  the  constant  of  the  Smagorinsky  model  is  adjusted  to  match  the  same  volume- 
averaged  production  as  the  “exact”  term,  the  value  of  Ck  varies  as  shown  in  Table  3.  The 
Smagorinsky  constant  can  also  be  allowed  to  vary  with  x3.  If  Ck  is  adjusted  at  each  x3  to 
give  the  simulation  results,  the  modified  Ck  varies  with  x3  as  shown  in  Figure  18.  A  value  of 
0.05  gives  best  agreement  in  the  center  region  and  the  value  of  Ck  varies  significantly  closer 
to  the  walls.  Ck  in  the  side  regions  increases  with  filter  width. 

<Pgc  >  is  plotted  in  Figure  19  for  the  simulation  data  (labeled  “exact”)  and  for  the 
two  models  mentioned  above.  <  Pr+c  >  for  the  Smagorinsky  model  was  calculated  at  the 
different  filter  widths  with  the  value  of  Ck  in  Table  3.  Since  the  scale  similarity  model  has 
no  adjustable  constant,  {Pr+c}  calculated  using  this  model  is  different  from  the  other  two 
in  the  horizontal  average  plots.  {-Pr+c)  versus  filter  width  is  plotted  in  Figure  20  where  the 
same  value  of  Ck  (0.07)  was  used  for  all  filter  widths.  The  data  in  Figure  19  gives  a  more 
local  comparison  of  the  models  to  the  “exact”.  At  the  small  filter  widths,  the  scale-similarity 
model  gives  good  agreement  with  the  “exact”  value  as  x3  (the  inhomogenous  direction)  is 
varied.  For  larger  filter  widths,  the  agreement  is  not  as  good  and  the  scale  similarity  model 
even  predicts  negative  values  for  <  Pr+c  >  in  the  core  region  while  the  “exact”  value 
remains  positive.  The  Smagorinsky  model  shows  no  agreement  with  the  variation  of  spacial 
distance  in  an  inhomogenous  direction  and  can  only  be  providing  a  global  agreement  in 
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filter  width  Smagorinsky  scale-similarity 


1/16 

1/8 

1/4 

1/2 

3/4 


0.07 

0.12 

0.19 

0.19 

0.14 


0.99 

0.97 

0.89 

0.74 

0.68 


Table  4:  Volume  average  of  the  correlation  of  Pr+c  calculated  directly  with  Pr+c  calculated 
using  a  model  for  Hij  and  Z7g. 


energy  dissipation  as  has  generally  been  believed.  For  comparison  of  these  results  with  other 
calculations,  the  fraction  of  energy  below  the  filter  width  shown  in  Table  2  may  be  a  better 
parameter  than  the  value  of  the  filter  width  itself. 

Previous  studies  have  estimated,  that  the  scale-similarity  model  has  a  very  small  contri¬ 
bution  to  the  subgrid  turbulent  production  and  dissipation  [1].  The  current  results  show 
that  the  fraction  of  the  “exact”  production  energy  predicted  by  the  scale-similarity  model 
varies  significantly  with  A /  (Figure  20).  For  the  smallest  A/  tested  (too  small  to  be  of  use 
in  a  LES),  the  scale-similarity  model  prediction  of  <  Pr+c  >  is  almost  equal  to  the  “exact” 
value  (Figure  19).  As  A /  increases,  this  fraction  predicted  decreases  to  12%  for  A/  =  0.75. 
As  mentioned  above,  the  value  of  Cr  that  best  matches  {Pr+c)  °f  the  Smagorinsky  model 
to  the  “exact”  value  also  varies  with  A/.  In  most  LES  studies  the  selection  of  a  model  and 
the  model  constants  is  reported  without  reference  to  the  fraction  of  energy  in  the  subgrid 
field.  Based  on  these  results  it  appears  that  the  choice  of  models  and  model  constants  should 
be  a  function  of  the  energy  fraction  in  the  subgrid  field. 

The  correlations  of  the  two  models  with  the  “exact”  term  based  on  the  simulation  data 
for  <  Ptf+C  >  are  shown  in  Figure  21.  The  horizontal  correlations  are  formed  as  follows: 


corr(x3)  = 


<  (a-  <  a  >)(&-  <b>)> 
yj<  (a-  <  a  >y  >J<  (<b-<b  >)a  > 


(65) 


where  a  and  b  are  functions  of  xi,  x%  and  x3.  The  correlations  of  Pr+c  and  Wr+c  (not 
shown)  have  roughly  the  same  magnitudes.  The  correlation  of  Wr+c  is  the  same  “scalar 
level”  comparison  in  the  study  done  by  Bardina  et.  al.[l],  except  that  the  correlations  were 
done  in  horizontal  planes  in  the  present  study.  The  vertical  average  of  the  results  in  Figure  21 
are  shown  in  Table  4.  The  correlation  levels  in  the  current  study  show  a  significantly  higher 
correlation  for  the  scale-similarity  model  than  for  the  Smagorinsky  model,  similar  to  that  in 
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[1].  In  [1]  the  correlation  coefficients  were  0.50  and  0.58  for  the  scale  similarity  model  and 
0.36  and  0.05  for  the  Smagorinsky  model.  The  flows  used  in  that  study  were  homogenous 
isotropic  turbulence  and  homogenous  shear  turbulence.  The  current  correlation  results  are 
in  general  agreement  with  [1].  The  scale  similarity  model  gives  better  local  agreement  with 
the  “exact"  results  than  the  Smagorinsky  model. 

8  Concluding  Remarks 

A  formal  study  of  the  set  of  kinetic  energy  equations  derived  from  filtering  the  original 
momentum  equations  has  provided  insight  into  the  physics  of  the  current  flow,  hi  particular, 
the  variation  of  the  thermal  and  velocity  production  terms  with  scale  (or  A/)  has  been 
quantified.  When  one  deals  with  the  large  quantity  of  detailed  data  that  exist  in  turbulent 
flows  some  averaging  of  the  details  is  necessary  to  study  the  flow.  Fourier  spectral  analysis 
provides  information  on  the  specially  averaged  energy  flow.  The  filter  analysis  presented  here 
provides  similar  information  but  with  a  specially  local  average.  Also,  the  current  study  has 
shown  that  the  diffusion  terms  are  not  negligible,  especially  near  the  walls.  The  Smagorinsky 
models  assumes  production  equals  dissipation,  neglecting  diffusion. 

The  extension  of  the  calculation  of  the  terms  of  the  filtered  set  of  kinetic  energy  equations 
to  evaluate  subgrid  models  is  natural.  One  difference  in  this  study  from  previous  “a  priori" 
studied  is  that  the  fine  scale  calculations  of  the  direct  simulation  were  not  projected  onto  a 
coarse  grid  before  calculating  the  model  results.  This  was  done  to  decouple  the  effect  of  the 
subgrid  filtering  from  the  filtering  that  results  from  projecting  onto  a  couse  grid. 

The  validity  of  the  “a  priori"  studies  to  investigate  LES  models  using  DNS  data  is  an 
open  issue.  Possibly  the  approach  taken  in  this  study  when  applied  to  a  wider  range  of  flows 
and  then  compared  in  actual  LES  studies  can  help  clarify  the  process  of  determining  subgrid 
models.  In  particular,  a  better  understanding  of  the  effect  of  combining  the  Smagorinsky 
model  with  the  scale  similarity  model  (called  the  linear  combination  model)  is  needed.  In 
the  study  where  the  linear  combination  model  was  proposed  [1],  it  was  suggested  that  the 
scale  similarity  part  contributed  little  (about  5%)  to  the  energy  dissipation.  The  current 
results  suggest  that  the  amount  varies  with  A/  and  is  generally  a  larger  percentage. 

9  Appendix  A  -  The  Rayleigh-Benard  problem 

The  Rayleigh-Benard  problem  is  a  simple  geometry,  laboratory-type  problem  used  to  study 
natural  convection.  The  problem  entails  the  study  of  the  fluid  motion  and  thermal  convection 
of  a  rectangular  fluid  layer  which  is  typically  heated  from  below.  The  layer  is  typically 
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thin,  has  no  mean  pressure  gradient  in  the  horisontal  directions,  and  has  uniform  boundary 
conditions  in  the  horisontal  directions. 

Wall  boundary  conditions  are  assumed  at  the  upper  and  lower  fluid  edge*.  Periodic 
boundary  conditions  are  assumed  in  the  horisontal  directions.  The  vertical  coordinate,  x3, 
was  defined  as  sero  halfway  between  the  two  walls  for  this  study. 

1.  Variables: 

•  fluid  properties  (all  constant)  - 

v  kinematic  viscosity 

a  thermal  diffusivity 

fi  coefficient  of  thermal  expansion 
p0  reference  fluid  density 
g  acceleration  of  gravity 

•  geometry  constants  - 

h  distance  between  the  two  horizontal  walls 

•  dependent  variables  - 

u i  fluid  velocity,  a  =  1,2, 3 
Tm  fluid  temperature 
Pm  fluid  pressure 
p  fluid  density 

•  independent  variables  * 

Xi  special  coordinates,  »  =  1,2,3 
t  time 

•  flow  constants  - 

Ta  temperature  at  lower  wall 

AT  temperature  difference  between  upper  and  lower  wall 

2.  Boundary  conditions: 

=  0  x3  =  ±h/2 

Ta  =  T0  x3  =  -h/2 

Ta  =  T0  -  AT  x3  =  h/2 
AT  >  0 

3.  Boussinesq  assumption: 

P9  =  Po9  -  po9P(Ta  -  T0)  (66) 

4.  Removal  of  static  temperature  and  pressure  gradient 

T.(2,t)  =  T,(z,t)  +  %  (67) 
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wo — wo — -££-(* $ + h) 

(88) 

dpm  dp  . 

a*i“  ax<+*?>  *"1,2 

(69) 

ap.  ap  /  at,  a 

“feT  — dT,+M{1+  + *>) 

(70) 

Non-dimeasionalise  by  a,  AT,  p0  and  h 

.  - 

114  ~o  /h 

(71) 

P  -  P 

Po<*IV 

(72) 

fs  — 

AT 

(73) 

va 

(74) 

Pr~  — 
a 

(75) 

*.«a 

(76) 

(77) 

Governing  Equations 

dtb 

—  =  0 
ax,- 

(78) 

(79) 

ar  ,  .  ar  a»f  A 
at  +ttia*,~a*7+1i3 

(80) 

fl  A  0* 

a  a&  duj 

<i""a^  +  a^ 

(81) 

In  the  main  body  of  this  paper  the  non-dimensional  variables  and  equations  will  be  used. 
The A  notation  will  be  removed  for  convenience. 


22 


References 


[1]  Bardina,  J.,  Ferziger,  J.H.  and  Reynolds,  W.C.,  1983,  Improved  turbulence  models 
based  on  large  eddy  simulation  of  homogeneous,  incompressible,  turbulent  flows,  Report 
No.  TF-19,  Dept,  of  Mech.  Engr.,  Standford  University. 

[2]  Carroll,  J.J.,  1976,  The  thermal  structure  of  turbulent  convection.  J.  Atmos.  Sc.  33, 
642. 

[3]  Busse,  F.H.,  1985,  Transition  to  turbulence  in  Rayleigh-  Benard  Convection.  Topics  in 
Applied  Physics:  Hydrodynamic  Instabilities  and  the  Transition  to  Turbulence  45,  ed. 
H.L.  Swinney  and  J.P.  Gollub,  Springer- Verlag. 

[4]  Clark,  R.A.,  Ferziger,  J.H.  and  Reynolds,  W.C.,  1979,  Evaluation  of  subgrid-scale  tur¬ 
bulence  models  using  a  fully  simulated  turbulent  flow.  J.  Fluid  Mech.  91,  1. 

[5]  Deardorff,  J.W.  and  Willis,  G.E.,  1967,  Investigation  of  turbulent  convection  between 
horizontal  plates.  J.  Fluid  Mech.  28,  675. 

[6]  Deardorff,  J.W.,  1973,  The  use  of  subgrid  transport  equations  in  a  three-dimensional 
model  of  atmospheric  turbulence.  J.  Fluid  Engineering  95,  429. 

[7]  Eidson,  T.M.,  1985,  Numerical  simulation  of  the  turbulent  Rayleigh-Benard  problem 
using  subgrid  modeling.  J.  Fluid  Mech.  158, 145. 

[8]  Eidson,  T.M.,  Hussaini,  M.Y.  and  Zang,  T.A.,  1986,  Simulation  of  the  turbulent 
Rayleigh-Benard  problem  using  a  spectal/finite  difference  technique.  In  Proceeding  of  the 
EUROMECH  Colloquium  No.  199  (ed..  U.  Schumann  and  R.Friedrich),  p.  188,  Vieweg, 
Braunschweig. 

[9]  Grotzbach,  G.,  1982,  Direct  numerical  simulation  of  laminar  and  turbulent  Benard 
convection.  J.  Fluid  Mech.  119,  145. 

[10]  Lilly,  D.K.,  1967,  The  representation  of  small-scale  turbulence  in  numerical  simulation 
experiments.  In  Proc.  of  the  IBM  Scientific  Computer  Symposium  on  Environmental 
Sciences ,  IBM  Form  No.  320-1951,  p.  195. 

[11]  Piomelli,  U.,  Zang,  T.A.,  Speziale,  C.G.  and  Hussaini,  M.Y.,  1989,  On  the  large-eddy 
simulation  of  transition  wall-bounded  flows.  Phy3.  Fluids  A  2,  257. 

[12]  Tennekes,  H.  and  Lumley,  J.L.,  1972,  A  First  Course  in  Turbulence ,  M.I.T.  Press, 
Cambridge,  Mass. 


24 


[13]  Reynolda,  W.C.,  1989,  The  Potential  and  Limitation!  of  Direct  and  Large  Eddy  Simu¬ 
lations.  Lecture  Note »  in  Phyeiece:  Whither  Turbulence?  Turbulence  at  the  Grostroadt 
357,  ed.  J.L.  Lumley,  Springer-  Verlag. 

[14]  Rogallo,  R.S.  and  Moin,  P.,  1984,  Numerical  simulation  of  turbulent  flows.  Ann.  Rev . 
FI  Mech.  18, 99. 


25 


Figure  3:  Volume-averaged  terms  in  the  set  of  filtered  energy  equations 
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Figure  4:  Volume-averaged  energy  flow  diagram  for  filter  width  0.25  .  All 
data  has  been  divided  by  10s. 
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Figure  5:  Volume-averaged  energy  flow  diagram  for  filter  width  0.50  .  All 
data  has  been  divided  by  10e. 


30 


Figure  6:  Horizontally-averaged  terms  of  the  total  kinetic  and  potential 
ergy  equation 
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Figure  10:  Horizontally-averaged  terms  in  the  filtered  kinetic  energy  equation 


(a)  Non-zero  principal  terms 


(b)  Components  of  velocity  production 

Figure  11:  Volume-averaged  terms  in  the  filtered  kinetic  energy  equation 
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Figure  13:  Horizontally-averaged  terms  in  the  cross  kinetic  energy  equation 


Figure  15:  Horizontally-averaged  terms  in  the  fluctuation  kinetic  energy 
equation 
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Figure  16:  Horizontally-averaged  components  of  the  velocity  production  term 
in  the  fluctuation  kinetic  energy  equation 
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Figure  17:  Volume-averaged  terms  in  the  fluctuation  kinetic  energy  equation 
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Figure  20:  Volume  average  of  Pr+c  '  simulation  and  2  models 
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1>.  Abatract 

A  filtering  analysis  of  a  turbulent  flow  has  been  developed  which  provides 
details  of  the  path  of  the  kinetic  energy  of  the  flow  from  its  creation  via  thermal 
production  to  its  dissipation.  A  low-pass  spatial  filter  is  used  to  split  the  velo¬ 
city  and  the  temperature  field  into  a  filtered  component  (composed  mainly  of  scales 
larger  than  a  specific  size,  nominally  the  filter  width)  and  a  fluctuation  component 
(scales  smaller  than  a  specific  size) .  Variables  derived  from  these  fields  can  fall 
into  one  of  the  above  two  ranges  or  be  composed  of  a  mixture  of  scales  dominated  by 
scales  near  the  specific  size.  The  filter  is  used  to  split  the  kinetic  energy  equa¬ 
tion  into  three  equations  corresponding  to  the  three  scale  ranges  described  above. 

The  data  from  a  direct  simulation  of  the  Rayleigh-Benard  problem  for  condi¬ 
tions  where  the  flow  is  turbulent  is  used'  to  calculate  the  individual  terms  in  the 
three  kinetic  energy  equations.  This  is  done  for  a  range  of  filter  widths.  These 
results  are  used  to  study  the  spatial  location  and  the  scale  range  of  the  thermal 
engergy  production,  the  cascading  of  kinetic  energy,  the  diffusion  of  kinetic  energy 
and  the  energy  dissipation.  These  results  are  used  also  to  evaluate  two  subgrid  mod¬ 
els  typically  used  in  large-eddy  simulations  of  turbulence.  Subgrid  models  attempt 
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